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ABSTRACT 



o 



Context. Atomic diffusion has been recognized as an important process that has to be considered in any computations of stellar 
models. In solar-type and cooler stars, this process is dominated by gravitational settling, which is now included in most stellar 
evolution codes. In hotter stars, radiative accelerations compete with gravity and become the dominant ingredient in the diffusion flux 
for most heavy elements. Introducing radiative accelerations into the computations of stellar models modifies the internal element 
distribution and may have major consequences on the stellar structure. Coupling these processes with hydrodynamical stellar motions 
has important consequences that need to be investigated in detail. 

Aims. We aim to include the computations of radiative accelerations in a stellar evolution code (here the TGEC code) using a simplified 
method (SVP) so that it may be coupled with sophisticated macroscopic motions. We also compare the results with those of the 
Montreal code in specific cases for validation and study the consequences of these coupled processes on accurate models of A- and 
I ■ early-type stars. 

Methods. We implemented radiative accelerations computations into the Toulouse-Geneva stellar evolution code following the semi- 
analytical prescription proposed by Alecian and LeBlanc. This allows more rapid computations than the full description used in the 
CZ2 . Montreal code. 

Results. We present results for A-type stellar models computed with this updated version of TGEC and compare them with similar 
published models obtained with the Montreal evolution code. We discuss the consequences for the coupling with macroscopic motions, 
including thermohaline convection. 
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chemical transport process. 

Accurate stellar modeling has recently been given a new boost Observations of chemical abundance anomalies in stars and 



(N ! with the advent of asteroseismology. The observations of os- unexpected stellar seismic behaviors have proved the necessity 

■ cillating stars and the analysis of the stellar oscillation proper- of ^eluding "non standard processes" in stellar evolutionary 

> i ties brought new and powerful constraints on stellar models and computations. In this framework, the main challenges encoun- 

• ■ allowed major progress in our understanding of stellar internal tered toda y b Y stellar physicists are to better determine the ef- 

X 1 structur e. Furthermore, sinc e the di scovery of the firs t extrasolar fects of rotatlon and magnetic fields and understand the chemi- 

H ! planets (IWolszczan & Fraill fl99l iMavor & Oueto3fl995l) . the cal transport processes better. This last point, and more specifi- 

P. ■ spectacular development of the exoplanet research field has also call y modeling of atomic diffusion including the radiative accel- 

sparked renewed interest in stellar physics, the accurate knowl- eratlons on individual elements, represents a key ingredient for 

edge of the host star being a necessary condition for characteriz- accurate stellar modeling. 

ing the surrounding planets. The importance of atomic diffusion inside stars is well es- 
tablished: notonly can it modify the atmospheric abundances 

One of the most important successes of astrophysics was the (Mic haud1ll970h . as observed in the so-called "chemically pecu- 

understanding of the basics of stellar internal structure and evo- liar stars", but it can als o have strong impli cations for the stel- 

lution. This progress was supported by the computations of nu- lar internal structure (e.g. Ric hard et al.120011) . In main-sequence 

merical models. Growing computational resources contributed solar-type and cooler stars (below about 1.2 M ), the radiative 

to refining the description of the physics of the stellar medium accelerations on the heavy elements are generally slower than 

(equation of state, opacities, nuclear reactions, etc) and allowed gravity in absolute v alue, owing to the s mall radiation flux com- 

building a simplified but efficient and widely used "standard pared to hotter stars (Michaud et alJl976l) . The elements heavier 

model". However, this standard model does not take the effects than hydrogen sink, even if the efficiency of this sinking may be 



of rotation and magnetic fields or the occurrence of accretion or modulated by the radiative accelerations dTurcotte et alJ 1998h . 
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One of the biggest success of helioseismology was to prove 
the importance of atomic diffusion in the Sun. Its introduc- 
tion in solar evolutionary models has significantly improved 
the agreement between the sound speed profile inside solar 
models and that deduced from helioseismic inversion tech- 
niques (e.g. Christensen-Dalsgaard et al. 1996; Richard et al. 
1996; iBrun et al l 119981; ITurcotte et all 119981; ISchlattll 120021 
(This agreement has however been spoiled by the new abun- 
dances proposed by Asplund et al. 2005, 2009.) Gravitational 
settling is now introduced in most stellar evolution codes. 

The effects of atomic diffusion on the seismic frequency 
of main-sequence, sola r-type star models have been studied 
by [Theado et al. (2005). They have shown that element segre- 
gation significantly alters the internal structure of the models 
and their oscillations frequencies. The frequency differences be- 
tween models with and without diffusion reach several micro- 
Hertz for stars with masses greater than 1.3M . 

In hotter stars, the radiative accelerations may become sig- 
nificantly stronger than gravity for many metals, which are 
pushed up. The variations with depth of the radiative acceler- 
ations of specific elements combined with the selective effects 
of gravitational settling can lead to element accumulation or de- 
pletion in various stellar layers. In A and B-type stars, iron-peak 
element accumulation appear in the Z-opacity bump (located at 
-200000 K). The induced opacity increase may l ead to local con- 
vection (IRicher et al.ll2000HRichard et al.ll200ll) . 

The iron-peak element accumulation in the opacity bump 
region can help trigger stellar pulsations, therefore improving 
the agreement between seismic observations and theoretical fre- 
quency spectra in many stars: e.g. in y Dor adus, Am, SPB, /3 
Cephei or sdB stars (see iTheado et al.ll2009l for a detailed dis- 
cussion of this subject). 

The effects of the radiative accelerations on the oscillation 
frequencies have been tested bv lEscobar et al.l(l2012l) . Very weak 
effects are observed for models with masses up to 1 .28M Q , but 
significant effects are expected for more massive stars. 

Progress in the understanding of the physics of early-type 
stars is limited by the radiative accelerations not being computed 
in most stellar evolutionary codes. Up to now, the Montreal stel- 
lar evolution code is the only one in which a complete, accurate, 
and consi stent treatment of radiative accelerations has been in- 
trodu ced (IRicher et al.ll2000l: iRichard et aDl200lb ITurcotte et all 
2000). The price to pay for this accuracy is that the heavy and 
CPU-time consuming computations of atomic processes do not 
allow additional sophisticated treatments of macroscopic mo- 
tions. In the Montreal code, turbulence is only treated as an ex- 
tremely simplified process with a turbulent diffusion coefficient 
proportional to a parameterized power of the density. 

In this context, we introduced into the Toulouse-Geneva 
Evolution Code (hereafter TGEC) new computations of the ra- 
diative accelerations on heavy ele ments, following the semi- 
analy tical prescription propose d by lAlecian & LeBlancl (2002) 
and iLeBlanc & Alecianl d2004l) . This prescription leads to fast 
but reasonably accurate computations, which represent a good 
compromise between accuracy and CPU-time consumption and 
allows coupling with macroscopic motions. Such a treatment of 
abundance variations inside early-type stars is necessary for a 
good understanding of these stars. 

In the present pape r, we do not introdu ce thermohaline con- 
vection as described in Thea dcTet alj J2009) because we want to 
present a detailed comparison with the Montreal code in which 
this specific physical process is not included. Very precise tests 
of the results obtained by the two codes for iron accumulation 
inside nearly identical models with similar physics are still un- 



derway. In particular, the computations of iron fluxes seem to 
lead to differences in some cases, which still have to be under- 
stood. These tests are beyond the scope of the present paper and 
will be presented elsewhere in the near future. Here we present a 
first step in the comparisons, namely the detailed study of the ra- 
diative accelerations obtained with the two methods. Some abun- 
dance profiles are only shown as indicators of the results that the 
TGEC code is presently able to obtain. 

In the following section, we explain in detail the major im- 
provements implemented in the TGEC code. In Sect. 3, we 
present a comparison between the Montreal and TGEC computa- 
tions for two similar stellar models. In Sect. 4, additional models 
are presented to illustrate the capabilities of the updated TGEC 
version and its application fields. Our conclusions are given in 
Sect. 5. 

2. New opacities and atomic diffusion 
computations in TGEC 

The T oulouse-Geneva stell ar evolution code is described in de- 
tail inlHui -Bon-Hoal d2008l) : however, the code has recently un- 
dergone major improvements that are reported in the following 
sections. The code can follow the time-dependent abundance 
variations of 21 species (12 elements and their main isotopes: 
H, 3 He, 4 He, 6 Li, 7 Li, 9 Be, 10 B, I2 C, I3 C, 14 N, I5 N, 16 0, 17 0, 
I8 0, 20 Ne, 22 Ne, 24 Mg, 25 Mg, 26 Mg, 40 Ca, and 56 Fe) in detail. 
The remaining metals are collected into an average species Z. 

2.1. Opacities 

In the TGEC code, the opacities are computed using the OPCD 
v3.3 codes and datcQ (Seaton 2005). They allow computating 
self-consistent Rosseland opacities taking the detailed composi- 
tion of the chemical mixture into account. The opacities are re- 
calculated by considering the abundance variations at each time 
step and at each mesh point. 

2.2. Atomic diffusion 

2.2.1. Computational methods 

The diffusion computations are based on the Boltzmann equa- 
tion for dilute collision-dominated plasma. At equilibrium, the 
solution of the equation is the Maxwellian distribution function. 
Transport properties in stars are computed considering small de- 
viations from the Maxwellian distribution. In this framework, 
two different formalism have been proposed to obtain approxi- 
mate solutions to the Boltzmann equations. 

The first method lies on the Chapman-Enskog theory 
(IChapman & Cowlindll97dl) . which assumes that the total dis- 
tribution function of a given species can be written as a conver- 
gent series, each term of the series representing successive ap- 
proximations to the distribution function. The formalism is first 
applied to a test element diffusing in the stellar plasma, taking 
the diffusion of electrons and the induced electric field into ac- 
count. The computations lead to a statistical "diffusion veloc- 
ity" of the test-element with respect to the main component of 
the plasma. In stars, this basic component is generally hydrogen, 
so that the diffusion of every element is computed with respect 
to neutral hydrogen and/or to protons. Meanwhile, the hydrogen 
abundance is renormalized to recover the abundance consistency 

1 The OPCD_3.3 packages is available on the following website: 
http://cdsweb.u-strasbg.fr/topbase/op.html 



2 



S. Theado et al.: The new Toulouse-Geneva Stellar Evolution Code including radiative accelerations of heavy elements 



and local stellar equilibrium. For the case of helium, which has a 
non-negligible abundance, corrections on the diffusion velocity 
as proposed by Mo ntmerle & Michaud ( 1 19761) are added. 

The second method has been developed by Burgers] (1 19691) . 
It is based on the Grad 13 moment approximation and the use 
of a Fokker-Planck collision term in the Boltzmann equation. In 
this approach, separate flow and heat equations for each com- 
ponent of a multicomponent mixture are solved simultaneously. 
This method provides a more convenient way for handling mul- 
ticomponent gases than the Chapman-Enskog one but is heavier 
to apply. In both methods, the diffusion coefficients can be writ- 
ten as functions of the collisions integrals that depend on the 
exact nature of the interaction between colliding particules. 

In TGEC, the atom i c diff usion is computed following 
the IChapman & Cowling! d 1970b formalism. We use the grav- 
itational and therm al diffusion coefficients as derived by 
Paqu ette et al.l dl986l) . who performed a detailed computation of 
the collision integrals for ionized elements diffusing in an ion- 
ized medium. In the case of the diffusion of neutral atoms in 
ions, or reverse, the p olarization of the neut rals is taken into ac- 
count as proposed by Micha ud et al.l d 19781) . For neutrals atoms 
diffusing in a neutral medium, the rigid sphere approximation is 
used. The average diffusing charge Z,- of the particules is com- 
puted by solving the Saha-Boltzmann equations. The partition 
functions are limited to the statistical weight of fundamental lev- 
els, and no correction for high density or degenerate electrons 
is introduced. The set of Saha-Boltzmann equations are solved 
by using a Newton-Raphson iterative scheme. When the com- 
putations of radiative accelerations are added to the diffusion 
computations, the Chapman-Enskog method is easier to handle 
than Burgers'. This situation helps treating cases with short time 
scales and rapid variations, as shown in Sect. 4. 



2.2.2. Radiative accelerations 

The radiative accelerations on C, N, O, Ne, Mg, Ca, and Fe have 
be en included in the TGEC co de following the improved version 
of iLeBlanc & Alecianl (120041) of the semi- analytical prescrip- 
tion proposed by lAlecian & LeBlanq d2002l) . This method al- 
lows very fast computation of radiative accelerations with a rea- 
sonable accuracy. Radiative accelerat ions due to boun d-bound 
dAlecianll 19851 lAlecian &"ArtruF l990) and bound-free dAlecianl 
[1991 " transitions are obtained using a parametric form of the ra- 
diative acceleration equation. The basic idea of this parametric 
method is to derive formula where the terms depending explic- 
itly on atomic data (such as gf values for instance) are separated 
from those depending on the stellar plasma and the abundances 
of the considered ion (with the aim of accounting for the satura- 
tion effects). In this framework, the radiative accelerations may 
be approximated by calculating a single value for each parameter 
found in the related equations. This i s the so-called single-valu ed 
parameter (or SVP) approximation (LeBla nc & Aleci an 2004). 

The interface of the SVP method with the models consists 
of a set of six parameters per ion, which allows to estimate the 
radiative acc eleration of each element through simple algebraic 
expressions ( LeBl anc & Aleci an 2004), for each time step of the 
run of the evolution code. These parameters are determined only 
at the beginning of the computation through interpolation as a 
function of the stellar mass, inside a pre-established grid. The 
computation of the total radiative acceleration for a given ele- 
ment with SVP also requires computing the ions' relative popu- 
lations. This is included in the set of the SVP numerical routines 
added to TGEC. 



The SVP approximation may be implemented in existing 
codes. There are much less data to process than for detailed ra- 
diative acceleration calculations because complete and detailed 
monochromatic opacities for each element are not needed. A 
new grid of S VP-parameters, well fitted to the stellar mass range 
considered i n this work, was computed following the procedure 
described in LeBlanc & Alecianl d2004l). An impleme ntation of 
the SVP method was first used in Thea do et al.l ((2009). 



3. Comparison with the Montreal computations 

In this section, we compare the TGEC results to those obtained 
with the Montreal code. For this purpose we compare the internal 
structures and the results of the diffusion computations for two 
models: one computed with the Montreal code, the other com- 
puted with TGEC. The Montreal model chosen fo r this compar- 
ison is the 5.3D50-3, 1.7M -model presented in Richar det al.l 
d2001l) . A similar model was computed with TGEC, including 
input physics and initial parameters as close as possible to the 
Montreal model. As a first step, we present comparisons of seven 
elements (listed in Table 1), which are especially important for 
A type stars. 

The detailed physi c s of the Montreal m odel can be found 
inlRichardet all d200ll) . iRicher et alJ (|2000), and iTurcotte et all 
d!998l) . The key ingredients are reported in the following sec- 
tions and compared to those introduced in the TGEC model. 



3.1. Input physics 

3.1 .1 . Basic input physics 

Here is a list that compares the basic input physics used in the 
TGEC and Montreal models: 



Equatio n of state: the Montreal model uses the CEFF equation 
of state (Chris tensen-Da lsgaard & Daeppen, 1992). The TGEC 
mode l is computed using the OPAL2001 ( Rogers & Navfonovl 
2002) equation. 



Opacities: in both models, the opacities are computed at every 
point in the star and at each evolution time step. The Montreal 
computations take the local abundance of 21 ch emical elements 
into account (see Table 1 o flTurcotte et al.ll9 98) using the OPAL 
monochromatic data. In TGEC, the opacities are computed as 
described in Sect. 12. ll considering the detailed abundance of the 
elements listed in Table Q] 

Nuclear reactions: t he Montreal code uses the nucl ear energy 
generation routine of Bahc all & Pinsonneaulll {T992). The nu- 
clear reaction rates used in TGE C are from the analytical formu- 
lae of the NACRE compilation dAngulolll999l) . 

Convection: in both Montreal and TGEC models, the delim- 
itation of the convective zones is based on the Schwarzschild 
criterion. The energy flux in the convection zones is computed 
following the Bohm-Vitense mixing length formalism. In the 
Montreal model, convective zone mixing is modeled as a diffu- 
sion process, described by a diffusion coefficient D m(x (cf. Sect. 
13.1.31 ). A mixing length approximation is used for D mtx , which 
always leads to high D ml - x values and very homogeneous convec- 
tive zones. In TGEC, convective zones are assumed to be instan- 
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Table 1. Initial chemical composition (in mass fraction) for 
models presented in Section 3. 



H 


0.70 


( 3 He+ 4 He) 


0.28 


12 C 


3.466 x 10 3 


14 N 


1.063 x 10 3 


16 Q 


9.645 x 10- 3 


40 Ca 


7.469 x 10- 5 


56 Fe 


1.436 x 10~ 3 



taneously homogenized. The mixing length parameter a is taken 
to be equal to 1.68 in both the TGEC and Montreal models. 

Initial chemical mixture: the initi al mixture used in the 
Montreal model is given in Table 1 of lTurcotte et al.l (Q998). The 
initial composition introduced in the TGEC model is given in 
Table Q] of this paper: the chemical elements followed in both 
models are introduced with the same initial abundance. 



3.1.2. Atomic diffusion 

In the Montreal code, a tom ic diffusion is comput ed as described 
in lRichard etaf] (|2001|) and lTurcotte et al.l (I1998I) . The Montreal 
calculations take the time-dependent variations of 28 species (in- 
cluding isotopes) into account: H, 3 He, 4 He, 6 Li, 7 Li, 9 Be, I0 B, 
n B, 12 C, I3 C, N, O, Ne, Na, Mg, Al, Si, P, S, CI, Ar, K, Ca, Ti, 
Cr, Mn, Fe, and Ni. Monochromatic OPAL data are used to cal- 
culate the Rosseland opacities and the radiative accelerations at 
each evolution time step. The abundances are updated at every 
iteration over the stellar structure. The diffusion coefficients and 
velocities are determined by solving the Burgers's flow equa- 
tions for ionized gases (lBurgersll969l) for all diffusing elements. 
The collisions integrals are from Paquette et al. ( 1986). 

As a first step, we include, in the TGEC model, the atomic 
diffusion (including radiative accelerations) of seven elements 
(eight species including isotopes), especially important for A 
type stars: H, 3 He, 4 He, 12 C, 14 N, ie O, 40 Ca, and 56 Fe. The 
atomic diffusion is computed as described in Sect. | 2.2| The diffu- 
sion v elocities are computed following fhe lChapman & Cowling! 
< fl970h formalism in the t est-atom approximati on, and the diffu- 
sion coefficients are from lPaquette et all (119861) . 

3.1 .3. Macroscopic transport of chemical elements 

Macroscopic transport in the Montreal code. 
In the Montreal code, macroscopic transport processes, includ- 
ing (semi)convection and turbulent mixing, are modeled as diffu- 
sion processes by adding a pure diffusion term to each element's 
diffusion velocity V pl : 



V i(macro) = -(D T + D mix )- 



dlnX; 



dr 



(1) 



where Dj and D ml x are turbulent diffusion coefficients with D m ; x 
representing the effects of convective and semiconvective mo- 
tions and Dj parametrizing turbulent transport. The variable X, 
represents the mass fraction of species i. 

In radiative layers D mlx = 0, while in convective zones it 
is computed using the mixing length theory (as stated in Sect. 
13.1.1b . In semiconvection zones, D m ; x is assumed proportional 
to (V - V a( /)/(Vi - V), where V and V fl j stands for the local and 



_10 
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o 

-5 



-4 

log (AM/M) 



-2 



Fig. 1. Turbulent diffusion coefficient in the 170Myr, 1.7M - 
model computed with TGEC. D tur t is constant from the surface 
down to log T =5.3. Just below this region, D, u ,b is by definition 
50 times larger than the He atomic diffusion coefficient (see Eq. 
|2). The vertical dashed line locates the region where log T = 5.3. 

adiabatic logarithmic temperature gradients and Vi is a func- 
tion of V„ = dJnu/dlnP. We refer the reader to Sect. 2.2 of 
Richa rd et al. for a more detailed description of the treat- 

ment of convection. 

In the 5.3D50-3-model, the turbulent transport is chosen 
strong enough to guarantee, during the whole evolution period 
considered, a complete mixing throughout the region between 
the surface and the point where log To = 5.3 (see Fig. 2 of 
Richa rd et al. I I2001I This mixing mimics the effects of a Fe 
convection zone except that it is imposed throughout evolution. 
Below log Tq = 5.3, the diffusion coefficient Dt o beys the fol- 
lowin g algebraic dependence on density (Eq. 1 . of Riche r et al.l 
2000): 



D T = coD(He )\ — 



(2) 



where po = p(To), u = 50, n=3 and 

D(He) = 3.3 x 10~ 15 r 2 - 5 /[4pln(l + 1.125 x 10~ 16 r 3 /p)]. (3) 

Macroscopic transport in the TGEC code. 
In the TGEC code, the macroscopic motions of a given species 
/ are parametrized, as in the Montreal code, by including a tur- 
bulent diffusion term to the diffusion velocity of the considered 
species: 



Vi(turb) = -D T 



dlnXi 
dr ' 



(4) 



The TGEC code does not include a D m , x term, since convection 
is treated as an instantaneous dilution. The diffusion coefficient 
Dt is chosen to mimic the mixing (Dj + D m j X ) introduced in 
the 5.3D50-3 Montreal model. Figure[TJrepresents this diffusion 
coefficient inside ou r model, it may be compared to Fig. 2 of 
Rich ard et al.l (12001 1) . From the surface down to log T = 5.3, Dj 
is chosen to homogenize the stellar material, below this region 
Dmrb is computed using the same expression as introduced in the 
Montreal code (i.e. following Eq. |2j. 

3.2. Results 

The Montreal and TGEC codes compute stellar evolution from 
pre-main sequence up to the subgiant branch. They were as- 
sumed to be homogeneous on the pre-main sequence, atomic dif- 
fusion, and turbulent transport were introduced at the beginning 
of the main-sequence phase. 
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Fig. 2. Evolutionary tracks of two 1.7M Q sequences computed 
with TGEC or with the Montreal code. T he dashed curve rep- 
resents the 5.3D50-3 model presented in Rich ard et al.l (l200lh 
(and computed with the Montreal code). The solid curve repre- 
sents a model computed with TGEC and including similar input 
physic. The black dots and crosses respectively represent the 30 
and 400Myr-models. 

Figure[2]compares the evolutionary tracks of the two models 
(for clarity the pre-main-sequence evolution is not shown). In 
spite of different treatments of some physical ingredients (e.g. 
equation of state, opacities, nuclear reactions, diffusion), the two 
tracks appear relatively close to each other. 

We first give a detailed comparison between the Montreal 
and TGEC calculations for the two 30Myr models. Figure [3] 
compares the internal structure of these two models. The relative 
differences in pressure, temperature, and density never exceed 
6% (and even less for P and p). The maximum relative differ- 
ences are observed near log(AM/M) ss -6.3, which corresponds 
to the base of the potentia l iron convective zone . 

Similarly to Fig. 7 of iRichard eTall d200ll) . Fig. H displays 
the radiative accelerations and the local gravity of the two mod- 
els. The discrepancy between gMTR and grcEC remains smaller 
than around 0.3 dex except for O in the upper layers (in the 
convection zone). We notice that the average accuracy of the 
SVP method is estimated to +0.1 dex dAlecian & LeBland2 002) 
with respect to the detailed computations of ISeatonl dl997l) . 
However, the SVP approximation is less accurate for light el- 
ements (lighter than O) than for heavy ones because this method 
uses parameters determined at the position of the maximum of 
each ion relative population. Because there are fewer ioniza- 
tion states for light elements, the gap in temperature between 
these maxima is larger for them. The comparison of Fig. [4] is 
fairly satisfactory, since radiative accelerations are computed in 
completely different ways and use different atomic and opacity 
databases. 

Figure[5]presents the diffusion velocities of several elements 
for the same models. Despite the different atomic diffusion pre- 
scriptions used, the diffusion velocities are quite close in most 



of the star. Figure [6] shows the abundance profiles in the two 30 
and 400Myr- models. At 30 Myr the profiles are quite close for 
all the species followed in TGEC. However, at 400 Myr discrep- 
ancies in the abundance profiles are observed. These differences 
are small for most elements (He, C, N, O, and Fe) except for 
Ca. As shown in Fig. [2] the positions of the two 400Myr mod- 
els in the HR diagram are slightly different, and are in particular 
more distant from each other than the 30Myr models. As a con- 
sequence, their internal structure is also expected to be different. 
In this context, it seems difficult to disentangle the discrepan- 
cies in the models due to structure variations from those caused 
specifically by the difference in diffusion calculations. 



4. Application field of the TGEC code 

In this section, we propose to demonstrate the abilities of the 
TGEC code and its application field. Models of A-F stars in- 
cluding atomic diffusion and minimal mixing (i.e. only a mild 
mixing below the convective zones to avoid steep and unrealis- 
tic composition gradients at the transition between radiative and 
conve ctive regions) have already been presented in Thead oet al.l 
(120091) . These models, which were evolved until the end of the 
main-sequence phase, have shown the ability of the code to com- 
pute complete evolutionary tracks including radiative levitation 
effects in the presence of minimal mixing. 

In this paper, we propose another test of the abilities of the 
code to manage rapid variations in the chemical composition 
and the thermal structure. For this purpose we present a set of 
models computed with TGEC with m asses ranging from 1.5 to 
2.5M . Like the models presented in ITheado et a l. (2009), they 
include atomic diffusion and a mild mixing below the convective 
zones. In the present models, the iron convective zone, when 
it appears, is assumed to be connected and completely mixed 
with the H/He convective envelope. These models were evolved 
from pre-main sequence up to hydrogen-core exhaustion. They 
were assumed to be chemically homogeneous on the pre-main 
sequence. Atomic diffusion was introduced at the beginning of 
the main-sequence phase. Because the present computations are 
done to test the TGEC code, we did not include here the thermo- 
halin e convection that m ust be added for comparison with real 
stars (ITheado et al.ll2009l) . 



4.1. Input physics 

As previous ly described, our models w ere computed using the 
OPAL2001 dRogers & Navfonovll2002l) equation of st ate. The 
nucle ar reactions were from the NACRE compilation (Angulo 
1999), and the opacities were computed as described in Sect. 
12.11 Atomic diffusion (including the radiative forces) was intro- 
duced as described in Sect. 12.21 for the following elements: H, 
3 He, 4 He, 12 C, I4 N, 16 O, 40 C a, and 56 Fe. The in itial metal mix- 
ture was the solar mixture o f lGrevesse & Noelsl d!993l) . The ini- 
tial mass fractions for the diffusing species are given in Table [3] 
Convection was computed using the mixing length theory with 
a mixing length parameter a — 1.8. The HI and Hell convective 
zones were assumed to be connected by overshooting and mixed 
together. As stated previously, the iron convective region that 
may appear during main-sequence evolution was also assumed 
connected and completely mixed with the surface convective re- 
gion. 

To avoid the appearance of sharp and nonphysical abun- 
dance gradients at the transition between radiative and convec- 
tive regions, we introduced mild mixing below the surface con- 
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Fig. 3. Comparison of the internal structures (pressure, temperature and density) of the two 1 .7M -models represented in Fig. [2] 
The comparison is shown at 30Myrs. MTR stands for the model computed with the Montreal stellar evolution code while TGEC 
stands for the model computed with the Toulouse Geneva stellar evolution code. 




OB 
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Fe 




log (AM/M) 



log (AM/M) Z 



Fig. 4. Radiative accelerations (in cgs unit) at 30Myr in the TGEC and Montrealmodels represented in Fig. [2] The solid curves rep- 
resents the accelerations in the TGEC model, while the dashed curve represent those in the Montreal model. The long-dashed curve 
represents the local gravity in the TGEC model. (The local gravity in the Montreal model is not represented as it is indistinguishable 
from that of the TGEC model.) 



vective zone. Thi s mixing was modeled as a diffusion process 
dSchatzmanll 19691) as described in Eq.|4]with 

D T = D czb exp (]n2^2* ) , (5) 

where D cz i, and r cz i, are the value of Dj and the value of the ra- 
dius at the base of the convective zone, respectively and A is the 



half width of the mixing region, D cz b and A are free parameters 
chosen to produce a mild mixing to a small extent. The value 
of D cz h is taken as equal to 2 x 10 cm 2 .s _1 , and A is fixed to 
0.5% of the radius below the surface convective region. Figure 
[7] presents this diffusion coefficient below the surface convec- 
tive zone of 1.7M models at two evolutionary stages. At 170 
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Fig. 5. Diffusion velocities (in cgs units) below the surface convective zone in the Montreal (dashed curves) and TGEC (solid curves) 
30Myr-models represented in Fig. [4] The vertical lines locate the bottom of the completely mixed region. The peaks observed for 
Ca and Fe correspond to sign changes in the diffusion velocities: the Ca diffusion velocity is positive for log(AM/M) < -6 and 
-3.5 < log(AM/M) < -2.5. The Fe diffusion velocity is positive for -7 < log(AM/M) < -4.5. 

Table 2. Characteristics of the models. 





Age (Myr) 


X c 


|AZ C |/Z C 




\AL\/L 


T eff (K) 


\AT eff \/T eff 




30.2 


0.6938 


0.03% 


8.289 


1.25% 


7975.66 


0.18% 


Mtgec 


30.2 


0.6940 




8.393 




7989.90 






400 


0.5778 


0.81 % 


8.951 


1.51% 


7659.49 


0.70% 


Mtgec 


400 


0.5731 




9.086 




7713.35 





Table 3. Initial chemical composition (in mass fraction) for 
models presented in Section 4. 



H 


0.7112 


( 3 He+ 4 He) 


0.2714 


12 C 


0.298 1447E-02 


l4 N 


0.9218931E-03 


16q 


0.8375229E-02 


40 Ca 


0.6283632E-04 


56p e 


0.1148674E-02 



Myr, the surface convective region includes the H and He con- 
vective zones (there is no Fe convective zone), and at 588 Myrs 
the external convective zone is composed of the H, He, and Fe 
convective zones assumed connected and mixed. 

The physics introduced in these models is close to that of the 
r30-3M, 1.5M -model presented in iRichard etail (l200ll) . Like 



our models, the r30-3M model includes mild mixing below the 
He convective zone. The turbulent diffusion coefficient is cho- 
sen as large as 10 4 cm.s _1 at the base of the convective zone. It 
then rapidly decreases with depth, varying like p -3 . Turbulent 
transport is chosen to be large enough to guarantee complete 
mixing between the Fe and He convective zones whenever the 
Fe convective zone appears. The initial metal mixture and the 
mixing length parameter are, however, different in the Montreal 
and TGEC models. For CPU-time consuming reasons, the evo- 
lution of the Montreal model was stopped at 89 Myrs, and the 
results for a 1.5M model are the only ones presented. TGEC 
models with masses up to 2.5M Q are computed along the whole 
main-sequence phase. Since the physics of our models and the 
r30-3M Montreal model are close but not similar, no detailed 
comparison between them can be carried out. However, qualita- 
tive similarities are underlined. 
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Fig. 6. Abundance profiles in the TGEC (solid curves) and Montreal (dashed curves) models shown in Fig. [2] as a function of the 
outer mass fraction. log(X/Xo) represents the ratio between the current element mass fraction and its initial values (at t=0). The blue 
and orange curves present the abundances at 30 Myr and 400 Myr after the ZAMS, respectively. 
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Fig. 7. Turbulent diffusion coefficient (according to Eq. [5]) in a 
1 .7M -model computed with TGEC at two evolutionary stages. 
The solid line represents the diffusion coefficient below the sur- 
face convective zone (due to the ionization of H and He) at 
170Myr, the dashed line represents the diffusion coefficient be- 
low the iron convective zone (assumed connected with the H and 
He ones) at 588Myrs, the vertical dotted line represents the base 
of the iron convective zone. 



4.2. Results 

Figure [8] displays the evolutionary tracks of our models. (For 
clarity pre-main-sequence evolution is not represented.) In the 
following, we describe, as a representative example, the results 



obtained for a 2.1M model. We focus on iron-diffusion related 
features. 

Figure |9]presents the profiles of various quantities inside the 
2. 1 M model along the main-sequence evolution. The left and 
middle columns illustrate early main-sequence evolution, the 
right hand column displays later evolutionary stages. 

On the main sequence, the competition between the gravi- 
tational settling and the radiative acceleration leads to iron ac- 
cumulations in the outer regions of the model. We note iron 
enrichment in the surface H/He convective zone and in the Z- 
bump region (at around log(AM/M) = -7 (T a 200000K). In 
this region, iron enrichment significantly increases the opacity, 
hence the radiative gradient, and when this radiative gradient ex- 
ceeds the adiabatic gradient, the so-called iron convective zone 
appears. Here it is assumed to be connected and mixed up to 
the H/He surface convective zone through overshooting. As a 
result, the bottom of the surface convective zone sinks abruptly 
to log(AM/M) =a -7. The dilution in this thick convective zone 
decreases the iron abundance in the opacity bump: the opacity 
and the radiative gradient consequently decrease, which leads to 
the disappearance of the iron convective zone. The surface con- 
vective zone then recedes. A new cycle starts as radiative forces 
proceed in accumulating iron in the opacity bump region. The 
middle column of Fig.|9]shows a second iron accumulation/deep 
convection phase. 

The alternation between convective and radiative episodes 
in the Z-bump proceeds from ^20 to ^ 1 15 Myr. A thick convec- 
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Fig. 9. Internal structure of the 2.1M -model represented in Fig. [8] The profiles are represented as a function of the outer mass 
fraction and at various evolutionary steps. Upper panel: iron abundance profiles, 56 Fe/ 56 Feo represents the ratio between the current 
Fe abundance and its initial value (cf. Table |3J. The various curves are defined in the upper panels and are valid for all the plots 
of the respective column. Middle panel: Rosseland opacity profiles. Lower panel: difference between the radiative and adiabatic 
gradients. 



tive zone then stays for the subsequent main-sequence evolution. 
The right hand column of Fig. [9] illustrates the internal structure 
variations of our model during the later evolutionary phases. 

Figure [10] shows the position of the bottom of the surface 
convective zone of our model during the main sequence. At the 
beginning of the main-sequence phase, the convective zone in- 
cludes the H and He ionization regions. After a few million years 
and because of the He gravitational settling, the Hell convective 
zone disappears leading to an abrupt recession of the convec- 
tive region (whose bottom subsequently lies at log(AM/M) 
-11 .25). After this first receding episode, the convective region 



undergoes a lOOMyr-period of rapid variations during which the 
convective depth oscillates between log(AM/M) —11.2 and 
log(AM/M) st -6.5. After this period, a thick surface convec- 
tive zone appears that slowly deepens inside the interior during 
the rest of the main-sequence phase. 

The lasting thick convective zone results from two features 
explained below. 

- Because of structural variations forced by the star's evo- 
lution, the iron accumulation needed to induce convection 
in the opacity bump varies along the main sequence. After 
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Fig. 10. Variation with time of the bottom of the surface convective zone in the 2.1M model represented in Fig. [8] Left panel: 
zoom on the 0-0. 12Gyr period. Right panel: zoom on the 0. 12-0. 8Gyr period. The dashed vertical lines are located respectively at 
25.9Myr, 32.6Myr, 1 17.1Myr, 254.3Myr, 380.3Myr, 576.0Myr, and 726.0Myr. 
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Fig. 8. Evolutionary tracks of 1.5 to 2.5M models computed 
with atomic diffusion (including radiative acceleration effects) 
and a mild mixing below the convective enveloppe. 



115 Myrs, in particular, smaller iron accumulations in the 
Z-bump are needed to initiate the iron convective zone. 
The succession of convective sinking and receding episodes, 
which occurs from 20 to ^115Myr, leads to rapid varia- 
tions in the surface iron abundance. As an illustration of this 
phenomenom, Fig. QT| displays the iron abundance evolu- 
tion over 30Myrs (including 5 convective sinking/receding 
episodes). The global effect of these episodes is illustrated in 
Fig- El which presents the time-dependent variations of the 
averaged iron surface abundance (averaged over 20 Myrs). 
The diffusion-induced iron enrichments (in the Z-bump but 
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Fig. 11. Time dependent variations in the Fe surface abun- 
dance in the 2.1M model presented in Fig. [8] The iron abun- 
dance is shown on a 30Myr period including 5 convective 
sinking/reced i ng ep isodes. This graph is similar to Fig. 16 of 
iRichard et alJ (l200ll) (see text for more details). 



also in the H/He convective zone) combined with the deep 
convective mixing episodes leads on average to an iron en- 
richment of the surface convective zone. For high values of 
this enrichment, the convective dilution in the presence of an 
iron convective zone may become unable to efficiently de- 
crease the iron abundance in the Z-bump and to suppress the 
iron convective zone. A thick convective zone then persists 
during the subsequent evolution. 
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Fig. 12. Iron surface abundance along main-sequence evolution 
for the 2. 1M model represented in Fig. [8] During the 0-120Myr 
period (i.e. when the external convective zone undergoes rapid 
and large variations), the iron abundance shown is an averaged 
value (over 20Myr). 



According to Fig. [12] the iron abundance in the surface con- 
vective zone continues to increase after 115 Myrs (because of 
radiative accelerations), while the convective zone slowly deep- 
ens. After 400Myr, the bottom of the convective zone reaches 
iron-underabundant regions, which causes as a slow decrease in 
the iron surface abundance. 

Figure[13]presents the convective zone extension (vs time) in 
all computed models. For an easy comparison with other masses, 
the 2.1 M model is shown in this figure. A succession of con- 
vective sinking and receding episodes is observed in all the com- 
puted models. This period of rapid oscillations of the convective 
depth is shorter for higher stellar masses. In the 1.5 M model, 
it persists during most of the main-sequence phase, whereas in 
the 2.5 M model, a thick convective zone appears in less than 
80 Myrs. 

The patterns observed during the early main-sequence phase 
of our mode l s are consistent with the results presented by 
iRichard et all d200lh for their r30-3M, 1.5M model. The alter- 
nation between radiative and convective episodes in the opac- 
ity bump that occurs in our models is similarly observed in the 
Montreal model as d escribed in Sect. 3 (p aragraph "Model with 
mixing episodes") of Richard ~et al] (l200l|). The comp arison be- 
tween our Fig. [TTI and Fig. 16 of IRichard et alj d2001l) (middle 
panel) also shows similar behaviour to the surface iron abun- 
dances. But while the Montreal code is limited by time consum- 
ing issues (their r30-3M, 1.5M model is evolved for less than 
100 Myr), TGEC is able to manage the rapid structural and com- 
positional variations due to the succession of convective sinking 
and receding episodes and allows computing complete main- 
sequence evolutionary tracks for masses up to at least 2.5M . 



5. Conclusion 

The implementation of radiative diffusion effects in stellar mod- 
eling turns out as a necessary condition for any computations of 
accurate models of F, A, and B type stars. Different ways of in- 
troducing these effects are considered by theoreticians, but few 
reliable methods are presently known. 

The most accurate method consists in computing the time- 
dependent variations of the chemical species present in the stel- 
lar mixture while using a precise description of the diffusion 
process. This includes an accurate determination of the radia- 
tive accelerations applied to each species and the resolution of 
the Burgers equations to obtain the diffusion velocities. Such a 
method has been successfully introduced in the Montreal code, 
which has been used for many years for the modeling of main- 
sequence F-A-B type stars. This accurate method allows precise 
determinations of abundance variations due to atomic diffusion, 
but it leads to highly time-consuming codes, which makes the 
computations of complete main-sequence evolutionary tracks 
difficult. 

The development of asteroseismology during the past 30 
years and, in particular, the discovery of a wide variety of pulsat- 
ing F, A, and B-type stars have brought to light the necessity for 
developing stellar modeling tools able to provide accurate mod- 
els adapted to seimic analysis for these stars. The seismic mod- 
eling of a target star is an iterative process that requires com- 
putating numerous models before getting the best one. In this 
context, using a code like the Montreal one appears particularly 
unadapted. 

This encouraged us to develop an alternative code that in- 
cludes a more flexible prescription of the radiative accelera- 
tions. In this context we included the radiative diffusion in 
the TGEC code using t h e sem i-anal ytical prescription proposed 
by lAlecian & LeBland (120021) and LeBlanc & Alecianl (2004) 
(called the SVP method). This method, less accurate than the one 
introduced in the Montreal code, has been proven to give good 
results with reasonable precision. With the most recent version 
of TGEC, this method has been successfully introduced for the 
first time in a stellar evolution code. In TGEC, the diffusion ve- 
locities are then computed following the Chapman & Cowling 
formalism in the test-atom approximation. The time-dependent 
variations in the chemical species are computed consistently. 

To validate our code, we did a comparison between a model 
computed with the Montreal code and a similar model computed 
with TGEC. (The two models include input physics that are as 
close as possible.) The comparison shows good agreement in 
the computations of the radiative acceleration for models with 
similar internal structures. 

Nickel is not presently included in the TGEC c omputations, 
as its atomic data are not included in TOPBase dCunto et al. 
| 1993[) . the atomic data base of the Opacity Project (iSeaton et al. 
119921) from which the atomic data used in the SVP approxima- 
tion are taken. This explains its absence in the SVP package. 
Its introduction in TOPBastQ as well as in the SVP package, is 
underway. 

To demonstrate both the abilities of the TGEC code and its 
flexibility as compared to the Montreal one, we presented sets of 
models that include atomic diffusion with small additional mix- 
ing. In this context, models including "m inimal mixing" were 
already presented in iTheado et al.1 (l2009h : in Sect. [4] we chose 
to present a new set of models (with masses ranging from 1.5 
to 2.5M ) including physical ingredients close to the r30-3M, 



2 Nickel is, however, included in the opacity calculation of Opacity 
Project in an approximate way. 
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Fig. 13. Variation in the surface convective zone extension along the main-sequence phase for the same models as presented in Fig. 
[8] The vertical axis represents the outer mass fraction, the solid curves represent the bottom of the external convective zone. 



1.5Mq model of iRichard et alJ d200lh . In these models, only 
mild mixing is introduced below the surface convective zone, 
and the iron convection region, when its appears, is assumed to 
be connected and mixed with the surface convective envelope. 
Under these assumptions, our models show similar behavior to 
the Montreal models, at least during the first lOOMyr (where the 
Montreal computations stop). In particular, the introduced physi- 
cal hypothesis leads to an alternation of convective and radiative 
episodes in the opacity bump region in both the Montreal and 
the TGEC models. We have also shown that the TGEC code al- 
lows computing complete evolutionary tracks (up to the subgiant 
branch) for A-F stars with masses up to 2.5M . It is able to man- 
age rapid variations in the chemical composition, under cover of 
mild mixing at the transition between convective and radiative 
regions. The effects of quasi-pure atomic diffusion can then be 
evaluated in the context of stellar evolutionary models. 

As discussed in the introduction, precise tests of the results 
obtained by the two codes for the case of iron accumulation in- 
side stars, when thermohaline convection is not taken into ac- 
count, still have to be performed, since differences in the maxi- 
mum iron value appear in some cases. Such tests are underway 
and will be presented elsewhere in the near future. 

The improvements brought to the ato mic diffusion aspects 
of the TGEC code presented here and in iTheado et al.l (2009) 
makes it an excellent test bed for stellar evolution and astero- 
seismology. Its application to various types of stars will certainly 
help us understand them better. 
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ABSTRACT 



o 



Context. Atomic diffusion has been recognized as an important process that has to be considered in any computations of stellar 
models. In solar-type and cooler stars, this process is dominated by gravitational settling, which is now included in most stellar 
evolution codes. In hotter stars, radiative accelerations compete with gravity and become the dominant ingredient in the diffusion flux 
for most heavy elements. Introducing radiative accelerations into the computations of stellar models modifies the internal element 
distribution and may have major consequences on the stellar structure. Coupling these processes with hydrodynamical stellar motions 
has important consequences that need to be investigated in detail. 

Aims. We aim to include the computations of radiative accelerations in a stellar evolution code (here the TGEC code) using a simplified 
method (SVP) so that it may be coupled with sophisticated macroscopic motions. We also compare the results with those of the 
Montreal code in specific cases for validation and study the consequences of these coupled processes on accurate models of A- and 
I ■ early-type stars. 

Methods. We implemented radiative accelerations computations into the Toulouse-Geneva stellar evolution code following the semi- 
analytical prescription proposed by Alecian and LeBlanc. This allows more rapid computations than the full description used in the 
CZ2 . Montreal code. 

Results. We present results for A-type stellar models computed with this updated version of TGEC and compare them with similar 
published models obtained with the Montreal evolution code. We discuss the consequences for the coupling with macroscopic motions, 
including thermohaline convection. 

^ ' Key words. Stars: evolution - Stars: chemically peculiar - Diffusion 

o 

■ 1. Introduction mass loss into account, and it considers convection as the only 



in 



o 



chemical transport process. 

Accurate stellar modeling has recently been given a new boost Observations of chemical abundance anomalies in stars and 



(N ! with the advent of asteroseismology. The observations of os- unexpected stellar seismic behaviors have proved the necessity 

■ cillating stars and the analysis of the stellar oscillation proper- of ^eluding "non standard processes" in stellar evolutionary 

> i ties brought new and powerful constraints on stellar models and computations. In this framework, the main challenges encoun- 

• ■ allowed major progress in our understanding of stellar internal tered toda y b Y stellar physicists are to better determine the ef- 

X 1 structur e. Furthermore, sinc e the di scovery of the firs t extrasolar fects of rotatlon and magnetic fields and understand the chemi- 

H ! planets (IWolszczan & Fraill fl99l iMavor & Oueto3fl995l) . the cal transport processes better. This last point, and more specifi- 

P. ■ spectacular development of the exoplanet research field has also call y modeling of atomic diffusion including the radiative accel- 

sparked renewed interest in stellar physics, the accurate knowl- eratlons on individual elements, represents a key ingredient for 

edge of the host star being a necessary condition for characteriz- accurate stellar modeling. 

ing the surrounding planets. The importance of atomic diffusion inside stars is well es- 
tablished: notonly can it modify the atmospheric abundances 

One of the most important successes of astrophysics was the (Mic haud1ll970h . as observed in the so-called "chemically pecu- 

understanding of the basics of stellar internal structure and evo- liar stars", but it can als o have strong impli cations for the stel- 

lution. This progress was supported by the computations of nu- lar internal structure (e.g. Ric hard et al.120011) . In main-sequence 

merical models. Growing computational resources contributed solar-type and cooler stars (below about 1.2 M ), the radiative 

to refining the description of the physics of the stellar medium accelerations on the heavy elements are generally slower than 

(equation of state, opacities, nuclear reactions, etc) and allowed gravity in absolute v alue, owing to the s mall radiation flux com- 

building a simplified but efficient and widely used "standard pared to hotter stars (Michaud et alJl976l) . The elements heavier 

model". However, this standard model does not take the effects than hydrogen sink, even if the efficiency of this sinking may be 



of rotation and magnetic fields or the occurrence of accretion or modulated by the radiative accelerations dTurcotte et alJ 1998h . 
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One of the biggest success of helioseismology was to prove 
the importance of atomic diffusion in the Sun. Its introduc- 
tion in solar evolutionary models has significantly improved 
the agreement between the sound speed profile inside solar 
models and that deduced from helioseismic inversion tech- 
niques (e.g. Christensen-Dalsgaard et al. 1996; Richard et al. 
1996; iBrun et al l 119981; ITurcotte et all 119981; ISchlattll 120021 
(This agreement has however been spoiled by the new abun- 
dances proposed by Asplund et al. 2005, 2009.) Gravitational 
settling is now introduced in most stellar evolution codes. 

The effects of atomic diffusion on the seismic frequency 
of main-sequence, sola r-type star models have been studied 
by [Theado et al. (2005). They have shown that element segre- 
gation significantly alters the internal structure of the models 
and their oscillations frequencies. The frequency differences be- 
tween models with and without diffusion reach several micro- 
Hertz for stars with masses greater than 1.3M . 

In hotter stars, the radiative accelerations may become sig- 
nificantly stronger than gravity for many metals, which are 
pushed up. The variations with depth of the radiative acceler- 
ations of specific elements combined with the selective effects 
of gravitational settling can lead to element accumulation or de- 
pletion in various stellar layers. In A and B-type stars, iron-peak 
element accumulation appear in the Z-opacity bump (located at 
-200000 K). The induced opacity increase may l ead to local con- 
vection (IRicher et al.ll2000HRichard et al.ll200ll) . 

The iron-peak element accumulation in the opacity bump 
region can help trigger stellar pulsations, therefore improving 
the agreement between seismic observations and theoretical fre- 
quency spectra in many stars: e.g. in y Dor adus, Am, SPB, /3 
Cephei or sdB stars (see iTheado et al.ll2009l for a detailed dis- 
cussion of this subject). 

The effects of the radiative accelerations on the oscillation 
frequencies have been tested bv lEscobar et al.l(l2012l) . Very weak 
effects are observed for models with masses up to 1 .28M Q , but 
significant effects are expected for more massive stars. 

Progress in the understanding of the physics of early-type 
stars is limited by the radiative accelerations not being computed 
in most stellar evolutionary codes. Up to now, the Montreal stel- 
lar evolution code is the only one in which a complete, accurate, 
and consi stent treatment of radiative accelerations has been in- 
trodu ced (IRicher et al.ll2000l: iRichard et aDl200lb ITurcotte et all 
2000). The price to pay for this accuracy is that the heavy and 
CPU-time consuming computations of atomic processes do not 
allow additional sophisticated treatments of macroscopic mo- 
tions. In the Montreal code, turbulence is only treated as an ex- 
tremely simplified process with a turbulent diffusion coefficient 
proportional to a parameterized power of the density. 

In this context, we introduced into the Toulouse-Geneva 
Evolution Code (hereafter TGEC) new computations of the ra- 
diative accelerations on heavy ele ments, following the semi- 
analy tical prescription propose d by lAlecian & LeBlancl (2002) 
and iLeBlanc & Alecianl d2004l) . This prescription leads to fast 
but reasonably accurate computations, which represent a good 
compromise between accuracy and CPU-time consumption and 
allows coupling with macroscopic motions. Such a treatment of 
abundance variations inside early-type stars is necessary for a 
good understanding of these stars. 

In the present pape r, we do not introdu ce thermohaline con- 
vection as described in Thea dcTet alj J2009) because we want to 
present a detailed comparison with the Montreal code in which 
this specific physical process is not included. Very precise tests 
of the results obtained by the two codes for iron accumulation 
inside nearly identical models with similar physics are still un- 



derway. In particular, the computations of iron fluxes seem to 
lead to differences in some cases, which still have to be under- 
stood. These tests are beyond the scope of the present paper and 
will be presented elsewhere in the near future. Here we present a 
first step in the comparisons, namely the detailed study of the ra- 
diative accelerations obtained with the two methods. Some abun- 
dance profiles are only shown as indicators of the results that the 
TGEC code is presently able to obtain. 

In the following section, we explain in detail the major im- 
provements implemented in the TGEC code. In Sect. 3, we 
present a comparison between the Montreal and TGEC computa- 
tions for two similar stellar models. In Sect. 4, additional models 
are presented to illustrate the capabilities of the updated TGEC 
version and its application fields. Our conclusions are given in 
Sect. 5. 

2. New opacities and atomic diffusion 
computations in TGEC 

The T oulouse-Geneva stell ar evolution code is described in de- 
tail inlHui -Bon-Hoal d2008l) : however, the code has recently un- 
dergone major improvements that are reported in the following 
sections. The code can follow the time-dependent abundance 
variations of 21 species (12 elements and their main isotopes: 
H, 3 He, 4 He, 6 Li, 7 Li, 9 Be, 10 B, I2 C, I3 C, 14 N, I5 N, 16 0, 17 0, 
I8 0, 20 Ne, 22 Ne, 24 Mg, 25 Mg, 26 Mg, 40 Ca, and 56 Fe) in detail. 
The remaining metals are collected into an average species Z. 

2.1. Opacities 

In the TGEC code, the opacities are computed using the OPCD 
v3.3 codes and datcQ (Seaton 2005). They allow computating 
self-consistent Rosseland opacities taking the detailed composi- 
tion of the chemical mixture into account. The opacities are re- 
calculated by considering the abundance variations at each time 
step and at each mesh point. 

2.2. Atomic diffusion 

2.2.1. Computational methods 

The diffusion computations are based on the Boltzmann equa- 
tion for dilute collision-dominated plasma. At equilibrium, the 
solution of the equation is the Maxwellian distribution function. 
Transport properties in stars are computed considering small de- 
viations from the Maxwellian distribution. In this framework, 
two different formalism have been proposed to obtain approxi- 
mate solutions to the Boltzmann equations. 

The first method lies on the Chapman-Enskog theory 
(IChapman & Cowlindll97dl) . which assumes that the total dis- 
tribution function of a given species can be written as a conver- 
gent series, each term of the series representing successive ap- 
proximations to the distribution function. The formalism is first 
applied to a test element diffusing in the stellar plasma, taking 
the diffusion of electrons and the induced electric field into ac- 
count. The computations lead to a statistical "diffusion veloc- 
ity" of the test-element with respect to the main component of 
the plasma. In stars, this basic component is generally hydrogen, 
so that the diffusion of every element is computed with respect 
to neutral hydrogen and/or to protons. Meanwhile, the hydrogen 
abundance is renormalized to recover the abundance consistency 

1 The OPCD_3.3 packages is available on the following website: 
http://cdsweb.u-strasbg.fr/topbase/op.html 
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and local stellar equilibrium. For the case of helium, which has a 
non-negligible abundance, corrections on the diffusion velocity 
as proposed by Mo ntmerle & Michaud ( 1 19761) are added. 

The second method has been developed by Burgers] (1 19691) . 
It is based on the Grad 13 moment approximation and the use 
of a Fokker-Planck collision term in the Boltzmann equation. In 
this approach, separate flow and heat equations for each com- 
ponent of a multicomponent mixture are solved simultaneously. 
This method provides a more convenient way for handling mul- 
ticomponent gases than the Chapman-Enskog one but is heavier 
to apply. In both methods, the diffusion coefficients can be writ- 
ten as functions of the collisions integrals that depend on the 
exact nature of the interaction between colliding particules. 

In TGEC, the atom i c diff usion is computed following 
the IChapman & Cowling! d 1970b formalism. We use the grav- 
itational and therm al diffusion coefficients as derived by 
Paqu ette et al.l dl986l) . who performed a detailed computation of 
the collision integrals for ionized elements diffusing in an ion- 
ized medium. In the case of the diffusion of neutral atoms in 
ions, or reverse, the p olarization of the neut rals is taken into ac- 
count as proposed by Micha ud et al.l d 19781) . For neutrals atoms 
diffusing in a neutral medium, the rigid sphere approximation is 
used. The average diffusing charge Z,- of the particules is com- 
puted by solving the Saha-Boltzmann equations. The partition 
functions are limited to the statistical weight of fundamental lev- 
els, and no correction for high density or degenerate electrons 
is introduced. The set of Saha-Boltzmann equations are solved 
by using a Newton-Raphson iterative scheme. When the com- 
putations of radiative accelerations are added to the diffusion 
computations, the Chapman-Enskog method is easier to handle 
than Burgers'. This situation helps treating cases with short time 
scales and rapid variations, as shown in Sect. 4. 



2.2.2. Radiative accelerations 

The radiative accelerations on C, N, O, Ne, Mg, Ca, and Fe have 
be en included in the TGEC co de following the improved version 
of iLeBlanc & Alecianl (120041) of the semi- analytical prescrip- 
tion proposed by lAlecian & LeBlanq d2002l) . This method al- 
lows very fast computation of radiative accelerations with a rea- 
sonable accuracy. Radiative accelerat ions due to boun d-bound 
dAlecianll 19851 lAlecian &"ArtruF l990) and bound-free dAlecianl 
[1991 " transitions are obtained using a parametric form of the ra- 
diative acceleration equation. The basic idea of this parametric 
method is to derive formula where the terms depending explic- 
itly on atomic data (such as gf values for instance) are separated 
from those depending on the stellar plasma and the abundances 
of the considered ion (with the aim of accounting for the satura- 
tion effects). In this framework, the radiative accelerations may 
be approximated by calculating a single value for each parameter 
found in the related equations. This i s the so-called single-valu ed 
parameter (or SVP) approximation (LeBla nc & Aleci an 2004). 

The interface of the SVP method with the models consists 
of a set of six parameters per ion, which allows to estimate the 
radiative acc eleration of each element through simple algebraic 
expressions ( LeBl anc & Aleci an 2004), for each time step of the 
run of the evolution code. These parameters are determined only 
at the beginning of the computation through interpolation as a 
function of the stellar mass, inside a pre-established grid. The 
computation of the total radiative acceleration for a given ele- 
ment with SVP also requires computing the ions' relative popu- 
lations. This is included in the set of the SVP numerical routines 
added to TGEC. 



The SVP approximation may be implemented in existing 
codes. There are much less data to process than for detailed ra- 
diative acceleration calculations because complete and detailed 
monochromatic opacities for each element are not needed. A 
new grid of S VP-parameters, well fitted to the stellar mass range 
considered i n this work, was computed following the procedure 
described in LeBlanc & Alecianl d2004l). An impleme ntation of 
the SVP method was first used in Thea do et al.l ((2009). 



3. Comparison with the Montreal computations 

In this section, we compare the TGEC results to those obtained 
with the Montreal code. For this purpose we compare the internal 
structures and the results of the diffusion computations for two 
models: one computed with the Montreal code, the other com- 
puted with TGEC. The Montreal model chosen fo r this compar- 
ison is the 5.3D50-3, 1.7M -model presented in Richar det al.l 
d2001l) . A similar model was computed with TGEC, including 
input physics and initial parameters as close as possible to the 
Montreal model. As a first step, we present comparisons of seven 
elements (listed in Table 1), which are especially important for 
A type stars. 

The detailed physi c s of the Montreal m odel can be found 
inlRichardet all d200ll) . iRicher et alJ (|2000), and iTurcotte et all 
d!998l) . The key ingredients are reported in the following sec- 
tions and compared to those introduced in the TGEC model. 



3.1. Input physics 

3.1 .1 . Basic input physics 

Here is a list that compares the basic input physics used in the 
TGEC and Montreal models: 



Equatio n of state: the Montreal model uses the CEFF equation 
of state (Chris tensen-Da lsgaard & Daeppen, 1992). The TGEC 
mode l is computed using the OPAL2001 ( Rogers & Navfonovl 
2002) equation. 



Opacities: in both models, the opacities are computed at every 
point in the star and at each evolution time step. The Montreal 
computations take the local abundance of 21 ch emical elements 
into account (see Table 1 o flTurcotte et al.ll9 98) using the OPAL 
monochromatic data. In TGEC, the opacities are computed as 
described in Sect. 12. ll considering the detailed abundance of the 
elements listed in Table Q] 

Nuclear reactions: t he Montreal code uses the nucl ear energy 
generation routine of Bahc all & Pinsonneaulll {T992). The nu- 
clear reaction rates used in TGE C are from the analytical formu- 
lae of the NACRE compilation dAngulolll999l) . 

Convection: in both Montreal and TGEC models, the delim- 
itation of the convective zones is based on the Schwarzschild 
criterion. The energy flux in the convection zones is computed 
following the Bohm-Vitense mixing length formalism. In the 
Montreal model, convective zone mixing is modeled as a diffu- 
sion process, described by a diffusion coefficient D m(x (cf. Sect. 
13.1.31 ). A mixing length approximation is used for D mtx , which 
always leads to high D ml - x values and very homogeneous convec- 
tive zones. In TGEC, convective zones are assumed to be instan- 
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Table 1. Initial chemical composition (in mass fraction) for 
models presented in Section 3. 



H 


0.70 


( 3 He+ 4 He) 


0.28 


12 C 


3.466 x 10 3 


14 N 


1.063 x 10 3 


16 Q 


9.645 x 10- 3 


40 Ca 


7.469 x 10- 5 


56 Fe 


1.436 x 10~ 3 



taneously homogenized. The mixing length parameter a is taken 
to be equal to 1.68 in both the TGEC and Montreal models. 

Initial chemical mixture: the initi al mixture used in the 
Montreal model is given in Table 1 of lTurcotte et al.l (Q998). The 
initial composition introduced in the TGEC model is given in 
Table Q] of this paper: the chemical elements followed in both 
models are introduced with the same initial abundance. 



3.1.2. Atomic diffusion 

In the Montreal code, a tom ic diffusion is comput ed as described 
in lRichard etaf] (|2001|) and lTurcotte et al.l (I1998I) . The Montreal 
calculations take the time-dependent variations of 28 species (in- 
cluding isotopes) into account: H, 3 He, 4 He, 6 Li, 7 Li, 9 Be, I0 B, 
n B, 12 C, I3 C, N, O, Ne, Na, Mg, Al, Si, P, S, CI, Ar, K, Ca, Ti, 
Cr, Mn, Fe, and Ni. Monochromatic OPAL data are used to cal- 
culate the Rosseland opacities and the radiative accelerations at 
each evolution time step. The abundances are updated at every 
iteration over the stellar structure. The diffusion coefficients and 
velocities are determined by solving the Burgers's flow equa- 
tions for ionized gases (lBurgersll969l) for all diffusing elements. 
The collisions integrals are from Paquette et al. ( 1986). 

As a first step, we include, in the TGEC model, the atomic 
diffusion (including radiative accelerations) of seven elements 
(eight species including isotopes), especially important for A 
type stars: H, 3 He, 4 He, 12 C, 14 N, ie O, 40 Ca, and 56 Fe. The 
atomic diffusion is computed as described in Sect. | 2.2| The diffu- 
sion v elocities are computed following fhe lChapman & Cowling! 
< fl970h formalism in the t est-atom approximati on, and the diffu- 
sion coefficients are from lPaquette et all (119861) . 

3.1 .3. Macroscopic transport of chemical elements 

Macroscopic transport in the Montreal code. 
In the Montreal code, macroscopic transport processes, includ- 
ing (semi)convection and turbulent mixing, are modeled as diffu- 
sion processes by adding a pure diffusion term to each element's 
diffusion velocity V pl : 



V i(macro) = -(D T + D mix )- 



dlnX; 



dr 



(1) 



where Dj and D ml x are turbulent diffusion coefficients with D m ; x 
representing the effects of convective and semiconvective mo- 
tions and Dj parametrizing turbulent transport. The variable X, 
represents the mass fraction of species i. 

In radiative layers D mlx = 0, while in convective zones it 
is computed using the mixing length theory (as stated in Sect. 
13.1.1b . In semiconvection zones, D m ; x is assumed proportional 
to (V - V a( /)/(Vi - V), where V and V fl j stands for the local and 
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Fig. 1. Turbulent diffusion coefficient in the 170Myr, 1.7M - 
model computed with TGEC. D tur t is constant from the surface 
down to log T =5.3. Just below this region, D, u ,b is by definition 
50 times larger than the He atomic diffusion coefficient (see Eq. 
|2). The vertical dashed line locates the region where log T = 5.3. 

adiabatic logarithmic temperature gradients and Vi is a func- 
tion of V„ = dJnu/dlnP. We refer the reader to Sect. 2.2 of 
Richa rd et al. for a more detailed description of the treat- 

ment of convection. 

In the 5.3D50-3-model, the turbulent transport is chosen 
strong enough to guarantee, during the whole evolution period 
considered, a complete mixing throughout the region between 
the surface and the point where log To = 5.3 (see Fig. 2 of 
Richa rd et al. I I2001I This mixing mimics the effects of a Fe 
convection zone except that it is imposed throughout evolution. 
Below log Tq = 5.3, the diffusion coefficient Dt o beys the fol- 
lowin g algebraic dependence on density (Eq. 1 . of Riche r et al.l 
2000): 



D T = coD(He )\ — 



(2) 



where po = p(To), u = 50, n=3 and 

D(He) = 3.3 x 10~ 15 r 2 - 5 /[4pln(l + 1.125 x 10~ 16 r 3 /p)]. (3) 

Macroscopic transport in the TGEC code. 
In the TGEC code, the macroscopic motions of a given species 
/ are parametrized, as in the Montreal code, by including a tur- 
bulent diffusion term to the diffusion velocity of the considered 
species: 



Vi(turb) = -D T 



dlnXi 
dr ' 



(4) 



The TGEC code does not include a D m , x term, since convection 
is treated as an instantaneous dilution. The diffusion coefficient 
Dt is chosen to mimic the mixing (Dj + D m j X ) introduced in 
the 5.3D50-3 Montreal model. Figure[TJrepresents this diffusion 
coefficient inside ou r model, it may be compared to Fig. 2 of 
Rich ard et al.l (12001 1) . From the surface down to log T = 5.3, Dj 
is chosen to homogenize the stellar material, below this region 
Dmrb is computed using the same expression as introduced in the 
Montreal code (i.e. following Eq. |2j. 

3.2. Results 

The Montreal and TGEC codes compute stellar evolution from 
pre-main sequence up to the subgiant branch. They were as- 
sumed to be homogeneous on the pre-main sequence, atomic dif- 
fusion, and turbulent transport were introduced at the beginning 
of the main-sequence phase. 
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Fig. 2. Evolutionary tracks of two 1.7M Q sequences computed 
with TGEC or with the Montreal code. T he dashed curve rep- 
resents the 5.3D50-3 model presented in Rich ard et al.l (l200lh 
(and computed with the Montreal code). The solid curve repre- 
sents a model computed with TGEC and including similar input 
physic. The black dots and crosses respectively represent the 30 
and 400Myr-models. 

Figure[2]compares the evolutionary tracks of the two models 
(for clarity the pre-main-sequence evolution is not shown). In 
spite of different treatments of some physical ingredients (e.g. 
equation of state, opacities, nuclear reactions, diffusion), the two 
tracks appear relatively close to each other. 

We first give a detailed comparison between the Montreal 
and TGEC calculations for the two 30Myr models. Figure [3] 
compares the internal structure of these two models. The relative 
differences in pressure, temperature, and density never exceed 
6% (and even less for P and p). The maximum relative differ- 
ences are observed near log(AM/M) ss -6.3, which corresponds 
to the base of the potentia l iron convective zone . 

Similarly to Fig. 7 of iRichard eTall d200ll) . Fig. H displays 
the radiative accelerations and the local gravity of the two mod- 
els. The discrepancy between gMTR and grcEC remains smaller 
than around 0.3 dex except for O in the upper layers (in the 
convection zone). We notice that the average accuracy of the 
SVP method is estimated to +0.1 dex dAlecian & LeBland2 002) 
with respect to the detailed computations of ISeatonl dl997l) . 
However, the SVP approximation is less accurate for light el- 
ements (lighter than O) than for heavy ones because this method 
uses parameters determined at the position of the maximum of 
each ion relative population. Because there are fewer ioniza- 
tion states for light elements, the gap in temperature between 
these maxima is larger for them. The comparison of Fig. [4] is 
fairly satisfactory, since radiative accelerations are computed in 
completely different ways and use different atomic and opacity 
databases. 

Figure[5]presents the diffusion velocities of several elements 
for the same models. Despite the different atomic diffusion pre- 
scriptions used, the diffusion velocities are quite close in most 



of the star. Figure [6] shows the abundance profiles in the two 30 
and 400Myr- models. At 30 Myr the profiles are quite close for 
all the species followed in TGEC. However, at 400 Myr discrep- 
ancies in the abundance profiles are observed. These differences 
are small for most elements (He, C, N, O, and Fe) except for 
Ca. As shown in Fig. [2] the positions of the two 400Myr mod- 
els in the HR diagram are slightly different, and are in particular 
more distant from each other than the 30Myr models. As a con- 
sequence, their internal structure is also expected to be different. 
In this context, it seems difficult to disentangle the discrepan- 
cies in the models due to structure variations from those caused 
specifically by the difference in diffusion calculations. 



4. Application field of the TGEC code 

In this section, we propose to demonstrate the abilities of the 
TGEC code and its application field. Models of A-F stars in- 
cluding atomic diffusion and minimal mixing (i.e. only a mild 
mixing below the convective zones to avoid steep and unrealis- 
tic composition gradients at the transition between radiative and 
conve ctive regions) have already been presented in Thead oet al.l 
(120091) . These models, which were evolved until the end of the 
main-sequence phase, have shown the ability of the code to com- 
pute complete evolutionary tracks including radiative levitation 
effects in the presence of minimal mixing. 

In this paper, we propose another test of the abilities of the 
code to manage rapid variations in the chemical composition 
and the thermal structure. For this purpose we present a set of 
models computed with TGEC with m asses ranging from 1.5 to 
2.5M . Like the models presented in ITheado et a l. (2009), they 
include atomic diffusion and a mild mixing below the convective 
zones. In the present models, the iron convective zone, when 
it appears, is assumed to be connected and completely mixed 
with the H/He convective envelope. These models were evolved 
from pre-main sequence up to hydrogen-core exhaustion. They 
were assumed to be chemically homogeneous on the pre-main 
sequence. Atomic diffusion was introduced at the beginning of 
the main-sequence phase. Because the present computations are 
done to test the TGEC code, we did not include here the thermo- 
halin e convection that m ust be added for comparison with real 
stars (ITheado et al.ll2009l) . 



4.1. Input physics 

As previous ly described, our models w ere computed using the 
OPAL2001 dRogers & Navfonovll2002l) equation of st ate. The 
nucle ar reactions were from the NACRE compilation (Angulo 
1999), and the opacities were computed as described in Sect. 
12.11 Atomic diffusion (including the radiative forces) was intro- 
duced as described in Sect. 12.21 for the following elements: H, 
3 He, 4 He, 12 C, I4 N, 16 O, 40 C a, and 56 Fe. The in itial metal mix- 
ture was the solar mixture o f lGrevesse & Noelsl d!993l) . The ini- 
tial mass fractions for the diffusing species are given in Table [3] 
Convection was computed using the mixing length theory with 
a mixing length parameter a — 1.8. The HI and Hell convective 
zones were assumed to be connected by overshooting and mixed 
together. As stated previously, the iron convective region that 
may appear during main-sequence evolution was also assumed 
connected and completely mixed with the surface convective re- 
gion. 

To avoid the appearance of sharp and nonphysical abun- 
dance gradients at the transition between radiative and convec- 
tive regions, we introduced mild mixing below the surface con- 
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Fig. 3. Comparison of the internal structures (pressure, temperature and density) of the two 1 .7M -models represented in Fig. [2] 
The comparison is shown at 30Myrs. MTR stands for the model computed with the Montreal stellar evolution code while TGEC 
stands for the model computed with the Toulouse Geneva stellar evolution code. 
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Fig. 4. Radiative accelerations (in cgs unit) at 30Myr in the TGEC and Montrealmodels represented in Fig. [2] The solid curves rep- 
resents the accelerations in the TGEC model, while the dashed curve represent those in the Montreal model. The long-dashed curve 
represents the local gravity in the TGEC model. (The local gravity in the Montreal model is not represented as it is indistinguishable 
from that of the TGEC model.) 



vective zone. Thi s mixing was modeled as a diffusion process 
dSchatzmanll 19691) as described in Eq.|4]with 

D T = D czb exp (]n2^2* ) , (5) 

where D cz i, and r cz i, are the value of Dj and the value of the ra- 
dius at the base of the convective zone, respectively and A is the 



half width of the mixing region, D cz b and A are free parameters 
chosen to produce a mild mixing to a small extent. The value 
of D cz h is taken as equal to 2 x 10 cm 2 .s _1 , and A is fixed to 
0.5% of the radius below the surface convective region. Figure 
[7] presents this diffusion coefficient below the surface convec- 
tive zone of 1.7M models at two evolutionary stages. At 170 
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Fig. 5. Diffusion velocities (in cgs units) below the surface convective zone in the Montreal (dashed curves) and TGEC (solid curves) 
30Myr-models represented in Fig. [4] The vertical lines locate the bottom of the completely mixed region. The peaks observed for 
Ca and Fe correspond to sign changes in the diffusion velocities: the Ca diffusion velocity is positive for log(AM/M) < -6 and 
-3.5 < log(AM/M) < -2.5. The Fe diffusion velocity is positive for -7 < log(AM/M) < -4.5. 

Table 2. Characteristics of the models. 





Age (Myr) 


X c 


|AZ C |/Z C 




\AL\/L 


T eff (K) 


\AT eff \/T eff 




30.2 


0.6938 


0.03% 


8.289 


1.25% 


7975.66 


0.18% 


Mtgec 


30.2 


0.6940 




8.393 




7989.90 






400 


0.5778 


0.81 % 


8.951 


1.51% 


7659.49 


0.70% 


Mtgec 


400 


0.5731 




9.086 




7713.35 





Table 3. Initial chemical composition (in mass fraction) for 
models presented in Section 4. 



H 


0.7112 


( 3 He+ 4 He) 


0.2714 


12 C 


0.298 1447E-02 


l4 N 


0.9218931E-03 


16q 


0.8375229E-02 


40 Ca 


0.6283632E-04 


56p e 


0.1148674E-02 



Myr, the surface convective region includes the H and He con- 
vective zones (there is no Fe convective zone), and at 588 Myrs 
the external convective zone is composed of the H, He, and Fe 
convective zones assumed connected and mixed. 

The physics introduced in these models is close to that of the 
r30-3M, 1.5M -model presented in iRichard etail (l200ll) . Like 



our models, the r30-3M model includes mild mixing below the 
He convective zone. The turbulent diffusion coefficient is cho- 
sen as large as 10 4 cm.s _1 at the base of the convective zone. It 
then rapidly decreases with depth, varying like p -3 . Turbulent 
transport is chosen to be large enough to guarantee complete 
mixing between the Fe and He convective zones whenever the 
Fe convective zone appears. The initial metal mixture and the 
mixing length parameter are, however, different in the Montreal 
and TGEC models. For CPU-time consuming reasons, the evo- 
lution of the Montreal model was stopped at 89 Myrs, and the 
results for a 1.5M model are the only ones presented. TGEC 
models with masses up to 2.5M Q are computed along the whole 
main-sequence phase. Since the physics of our models and the 
r30-3M Montreal model are close but not similar, no detailed 
comparison between them can be carried out. However, qualita- 
tive similarities are underlined. 
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Fig. 6. Abundance profiles in the TGEC (solid curves) and Montreal (dashed curves) models shown in Fig. [2] as a function of the 
outer mass fraction. log(X/Xo) represents the ratio between the current element mass fraction and its initial values (at t=0). The blue 
and orange curves present the abundances at 30 Myr and 400 Myr after the ZAMS, respectively. 
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Fig. 7. Turbulent diffusion coefficient (according to Eq. [5]) in a 
1 .7M -model computed with TGEC at two evolutionary stages. 
The solid line represents the diffusion coefficient below the sur- 
face convective zone (due to the ionization of H and He) at 
170Myr, the dashed line represents the diffusion coefficient be- 
low the iron convective zone (assumed connected with the H and 
He ones) at 588Myrs, the vertical dotted line represents the base 
of the iron convective zone. 



4.2. Results 

Figure [8] displays the evolutionary tracks of our models. (For 
clarity pre-main-sequence evolution is not represented.) In the 
following, we describe, as a representative example, the results 



obtained for a 2.1M model. We focus on iron-diffusion related 
features. 

Figure |9]presents the profiles of various quantities inside the 
2. 1 M model along the main-sequence evolution. The left and 
middle columns illustrate early main-sequence evolution, the 
right hand column displays later evolutionary stages. 

On the main sequence, the competition between the gravi- 
tational settling and the radiative acceleration leads to iron ac- 
cumulations in the outer regions of the model. We note iron 
enrichment in the surface H/He convective zone and in the Z- 
bump region (at around log(AM/M) = -7 (T a 200000K). In 
this region, iron enrichment significantly increases the opacity, 
hence the radiative gradient, and when this radiative gradient ex- 
ceeds the adiabatic gradient, the so-called iron convective zone 
appears. Here it is assumed to be connected and mixed up to 
the H/He surface convective zone through overshooting. As a 
result, the bottom of the surface convective zone sinks abruptly 
to log(AM/M) =a -7. The dilution in this thick convective zone 
decreases the iron abundance in the opacity bump: the opacity 
and the radiative gradient consequently decrease, which leads to 
the disappearance of the iron convective zone. The surface con- 
vective zone then recedes. A new cycle starts as radiative forces 
proceed in accumulating iron in the opacity bump region. The 
middle column of Fig.|9]shows a second iron accumulation/deep 
convection phase. 

The alternation between convective and radiative episodes 
in the Z-bump proceeds from ^20 to ^ 1 15 Myr. A thick convec- 
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Fig. 9. Internal structure of the 2.1M -model represented in Fig. [8] The profiles are represented as a function of the outer mass 
fraction and at various evolutionary steps. Upper panel: iron abundance profiles, 56 Fe/ 56 Feo represents the ratio between the current 
Fe abundance and its initial value (cf. Table |3J. The various curves are defined in the upper panels and are valid for all the plots 
of the respective column. Middle panel: Rosseland opacity profiles. Lower panel: difference between the radiative and adiabatic 
gradients. 



tive zone then stays for the subsequent main-sequence evolution. 
The right hand column of Fig. [9] illustrates the internal structure 
variations of our model during the later evolutionary phases. 

Figure [10] shows the position of the bottom of the surface 
convective zone of our model during the main sequence. At the 
beginning of the main-sequence phase, the convective zone in- 
cludes the H and He ionization regions. After a few million years 
and because of the He gravitational settling, the Hell convective 
zone disappears leading to an abrupt recession of the convec- 
tive region (whose bottom subsequently lies at log(AM/M) 
-11 .25). After this first receding episode, the convective region 



undergoes a lOOMyr-period of rapid variations during which the 
convective depth oscillates between log(AM/M) —11.2 and 
log(AM/M) st -6.5. After this period, a thick surface convec- 
tive zone appears that slowly deepens inside the interior during 
the rest of the main-sequence phase. 

The lasting thick convective zone results from two features 
explained below. 

- Because of structural variations forced by the star's evo- 
lution, the iron accumulation needed to induce convection 
in the opacity bump varies along the main sequence. After 
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Fig. 10. Variation with time of the bottom of the surface convective zone in the 2.1M model represented in Fig. [8] Left panel: 
zoom on the 0-0. 12Gyr period. Right panel: zoom on the 0. 12-0. 8Gyr period. The dashed vertical lines are located respectively at 
25.9Myr, 32.6Myr, 1 17.1Myr, 254.3Myr, 380.3Myr, 576.0Myr, and 726.0Myr. 
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Fig. 8. Evolutionary tracks of 1.5 to 2.5M models computed 
with atomic diffusion (including radiative acceleration effects) 
and a mild mixing below the convective enveloppe. 



115 Myrs, in particular, smaller iron accumulations in the 
Z-bump are needed to initiate the iron convective zone. 
The succession of convective sinking and receding episodes, 
which occurs from 20 to ^115Myr, leads to rapid varia- 
tions in the surface iron abundance. As an illustration of this 
phenomenom, Fig. QT| displays the iron abundance evolu- 
tion over 30Myrs (including 5 convective sinking/receding 
episodes). The global effect of these episodes is illustrated in 
Fig- El which presents the time-dependent variations of the 
averaged iron surface abundance (averaged over 20 Myrs). 
The diffusion-induced iron enrichments (in the Z-bump but 
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Fig. 11. Time dependent variations in the Fe surface abun- 
dance in the 2.1M model presented in Fig. [8] The iron abun- 
dance is shown on a 30Myr period including 5 convective 
sinking/reced i ng ep isodes. This graph is similar to Fig. 16 of 
iRichard et alJ (l200ll) (see text for more details). 



also in the H/He convective zone) combined with the deep 
convective mixing episodes leads on average to an iron en- 
richment of the surface convective zone. For high values of 
this enrichment, the convective dilution in the presence of an 
iron convective zone may become unable to efficiently de- 
crease the iron abundance in the Z-bump and to suppress the 
iron convective zone. A thick convective zone then persists 
during the subsequent evolution. 
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Fig. 12. Iron surface abundance along main-sequence evolution 
for the 2. 1M model represented in Fig. [8] During the 0-120Myr 
period (i.e. when the external convective zone undergoes rapid 
and large variations), the iron abundance shown is an averaged 
value (over 20Myr). 



According to Fig. [12] the iron abundance in the surface con- 
vective zone continues to increase after 115 Myrs (because of 
radiative accelerations), while the convective zone slowly deep- 
ens. After 400Myr, the bottom of the convective zone reaches 
iron-underabundant regions, which causes as a slow decrease in 
the iron surface abundance. 

Figure[13]presents the convective zone extension (vs time) in 
all computed models. For an easy comparison with other masses, 
the 2.1 M model is shown in this figure. A succession of con- 
vective sinking and receding episodes is observed in all the com- 
puted models. This period of rapid oscillations of the convective 
depth is shorter for higher stellar masses. In the 1.5 M model, 
it persists during most of the main-sequence phase, whereas in 
the 2.5 M model, a thick convective zone appears in less than 
80 Myrs. 

The patterns observed during the early main-sequence phase 
of our mode l s are consistent with the results presented by 
iRichard et all d200lh for their r30-3M, 1.5M model. The alter- 
nation between radiative and convective episodes in the opac- 
ity bump that occurs in our models is similarly observed in the 
Montreal model as d escribed in Sect. 3 (p aragraph "Model with 
mixing episodes") of Richard ~et al] (l200l|). The comp arison be- 
tween our Fig. [TTI and Fig. 16 of IRichard et alj d2001l) (middle 
panel) also shows similar behaviour to the surface iron abun- 
dances. But while the Montreal code is limited by time consum- 
ing issues (their r30-3M, 1.5M model is evolved for less than 
100 Myr), TGEC is able to manage the rapid structural and com- 
positional variations due to the succession of convective sinking 
and receding episodes and allows computing complete main- 
sequence evolutionary tracks for masses up to at least 2.5M . 



5. Conclusion 

The implementation of radiative diffusion effects in stellar mod- 
eling turns out as a necessary condition for any computations of 
accurate models of F, A, and B type stars. Different ways of in- 
troducing these effects are considered by theoreticians, but few 
reliable methods are presently known. 

The most accurate method consists in computing the time- 
dependent variations of the chemical species present in the stel- 
lar mixture while using a precise description of the diffusion 
process. This includes an accurate determination of the radia- 
tive accelerations applied to each species and the resolution of 
the Burgers equations to obtain the diffusion velocities. Such a 
method has been successfully introduced in the Montreal code, 
which has been used for many years for the modeling of main- 
sequence F-A-B type stars. This accurate method allows precise 
determinations of abundance variations due to atomic diffusion, 
but it leads to highly time-consuming codes, which makes the 
computations of complete main-sequence evolutionary tracks 
difficult. 

The development of asteroseismology during the past 30 
years and, in particular, the discovery of a wide variety of pulsat- 
ing F, A, and B-type stars have brought to light the necessity for 
developing stellar modeling tools able to provide accurate mod- 
els adapted to seimic analysis for these stars. The seismic mod- 
eling of a target star is an iterative process that requires com- 
putating numerous models before getting the best one. In this 
context, using a code like the Montreal one appears particularly 
unadapted. 

This encouraged us to develop an alternative code that in- 
cludes a more flexible prescription of the radiative accelera- 
tions. In this context we included the radiative diffusion in 
the TGEC code using t h e sem i-anal ytical prescription proposed 
by lAlecian & LeBland (120021) and LeBlanc & Alecianl (2004) 
(called the SVP method). This method, less accurate than the one 
introduced in the Montreal code, has been proven to give good 
results with reasonable precision. With the most recent version 
of TGEC, this method has been successfully introduced for the 
first time in a stellar evolution code. In TGEC, the diffusion ve- 
locities are then computed following the Chapman & Cowling 
formalism in the test-atom approximation. The time-dependent 
variations in the chemical species are computed consistently. 

To validate our code, we did a comparison between a model 
computed with the Montreal code and a similar model computed 
with TGEC. (The two models include input physics that are as 
close as possible.) The comparison shows good agreement in 
the computations of the radiative acceleration for models with 
similar internal structures. 

Nickel is not presently included in the TGEC c omputations, 
as its atomic data are not included in TOPBase dCunto et al. 
| 1993[) . the atomic data base of the Opacity Project (iSeaton et al. 
119921) from which the atomic data used in the SVP approxima- 
tion are taken. This explains its absence in the SVP package. 
Its introduction in TOPBastQ as well as in the SVP package, is 
underway. 

To demonstrate both the abilities of the TGEC code and its 
flexibility as compared to the Montreal one, we presented sets of 
models that include atomic diffusion with small additional mix- 
ing. In this context, models including "m inimal mixing" were 
already presented in iTheado et al.1 (l2009h : in Sect. [4] we chose 
to present a new set of models (with masses ranging from 1.5 
to 2.5M ) including physical ingredients close to the r30-3M, 



2 Nickel is, however, included in the opacity calculation of Opacity 
Project in an approximate way. 
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Fig. 13. Variation in the surface convective zone extension along the main-sequence phase for the same models as presented in Fig. 
[8] The vertical axis represents the outer mass fraction, the solid curves represent the bottom of the external convective zone. 



1.5Mq model of iRichard et alJ d200lh . In these models, only 
mild mixing is introduced below the surface convective zone, 
and the iron convection region, when its appears, is assumed to 
be connected and mixed with the surface convective envelope. 
Under these assumptions, our models show similar behavior to 
the Montreal models, at least during the first lOOMyr (where the 
Montreal computations stop). In particular, the introduced physi- 
cal hypothesis leads to an alternation of convective and radiative 
episodes in the opacity bump region in both the Montreal and 
the TGEC models. We have also shown that the TGEC code al- 
lows computing complete evolutionary tracks (up to the subgiant 
branch) for A-F stars with masses up to 2.5M . It is able to man- 
age rapid variations in the chemical composition, under cover of 
mild mixing at the transition between convective and radiative 
regions. The effects of quasi-pure atomic diffusion can then be 
evaluated in the context of stellar evolutionary models. 

As discussed in the introduction, precise tests of the results 
obtained by the two codes for the case of iron accumulation in- 
side stars, when thermohaline convection is not taken into ac- 
count, still have to be performed, since differences in the maxi- 
mum iron value appear in some cases. Such tests are underway 
and will be presented elsewhere in the near future. 

The improvements brought to the ato mic diffusion aspects 
of the TGEC code presented here and in iTheado et al.l (2009) 
makes it an excellent test bed for stellar evolution and astero- 
seismology. Its application to various types of stars will certainly 
help us understand them better. 
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